﻿/**************************************************************
 * 
 * Copyright © WMSigner Team <www.wmsigner.com> 2007-2008
 * Copyright © Dmitry Kukushkin <codecity@gmail.com> 2007-2008
 * All Rights Reserved
 * 
 **************************************************************/

// Литература:
// [1] Handbook of Applied Cryptography, by A. Menezes, P. van Oorschot, and S. Vanstone, CRC Press, 1996.
// [2] М. Вельшенбах, Криптография на C и C++ в действии. М.: Издательство Триумф, 2004.

// Примечание:
// Данная реализация — 100% имплементация алгоритма описанного в [1]. Для удобства сохранены названия переменных.

// Благодарности:
// 1. CreatorCray, который помог проекту, предоставив свою реализацию MontgomeryExponentiation на языке C++.
// 2. Авторам класса BigInteger: Ben Maurer, Chew Keong, Sebastien Pouliot, Pieter Philippaerts. Был рассмотрен и отклонен избранный ими метод нахождения степени по модулю (Montgomery Reduction).

using System;

namespace WebMoney.Cryptography
{
    internal static class Montgomery
    {
        private const int BITS = 32;
        private const uint BIT_MASK = 0x80000000; // битовая маска

        public static uint Inverse(uint value)
        {
            uint temp = (((value + 2) & 4) << 1) + value;

            temp *= 2 - value * temp;
            temp *= 2 - value * temp;
            temp *= 2 - value * temp;

            return (uint)-temp;
        }

        // [1] 14.36 Algorithm Montgomery multiplication
        // INPUT: integers
        //      m = (m[n-1] ... m[1] m[0]){b},
        //      x = (x[n-1] ... x[1] x[0]){b},
        //      y = (y[n-1] ... y[1] y[0]){b}
        //           with 0 <= x, y < m,
        //           R = b^n with gcd(m, b) = 1,
        //           and mQ = -m^1 mod b.
        // OUTPUT: x * y * R^-1 mod m.
        public static uint[] MontgomeryMultiplication(uint[] x, uint[] y, uint[] m, uint mQ)
        {
            int n = ArithmeticKernel.Significance(m);

            // 1. A = 0. (Notation: A = (a[n] a[n-1] ... a[1] a[0]){b})
            uint[] a = new uint[n + 1];

            // 2. For i from 0 to (n - 1) do the following:
            for (int i = 0; i < n; i++)
            {
                // 2.1 u_i = (a[0] + x[i] * y[0]) * mQ mod b.
                uint u = (uint)((a[0] + (ulong)x[i] * y[0]) * mQ);

                // 2.2 A = (A + x[i] * y + u_i * m) / b.
                ulong xy = (ulong)x[i] * (ulong)y[0];
                ulong um = (ulong)u * (ulong)m[0];

                ulong temp = (ulong)a[0] + (ulong)((uint)xy) + (uint)um;
                ulong carry = (xy >> BITS) + (um >> BITS) + (temp >> BITS);

                for (int pos = 1; pos < n; pos++)
                {
                    xy = (ulong)x[i] * (ulong)y[pos];
                    um = (ulong)u * (ulong)m[pos];

                    temp = (ulong)a[pos] + (ulong)((uint)xy) + (uint)um + (uint)carry;
                    carry = (xy >> BITS) + (um >> BITS) + (temp >> BITS) + (carry >> BITS);

                    a[pos - 1] = (uint)temp;
                }

                carry += a[n];

                a[n - 1] = (uint)carry;
                a[n] = (uint)(carry >> BITS);
            }

            // 3. If A >= m then A = A - m
            if (ArithmeticKernel.Compare(a, m) >= 0)
                ArithmeticKernel.Sub(ref a, m);

            // 4. Return (A).
            return a;
        }

        public static ushort[] MontgomeryExponentiationUint16(uint[] xValue, uint[] e, uint[] m)
        {
            uint[] uint32result = MontgomeryExponentiation(xValue, e, m);

            ushort[] uint16Result = new ushort[uint32result.Length * 2];

            Buffer.BlockCopy(uint32result, 0, uint16Result, 0, uint32result.Length * 4);

            ArithmeticKernel.Normalize(ref uint16Result);

            return uint16Result;
        }

        // [1] 14.94 Algorithm Montgomery exponentiation
        // INPUT:
        //      m = (m[l-1] ... m[0]){b},
        //      R = b^l,
        //      mQ = -m^1 mod b,
        //      e = (e[t] ... e[0]){2}
        //           with e[t] = 1,
        //           and an integer x, 1 <= x < m.
        // OUTPUT: x^e mod m.
        public static uint[] MontgomeryExponentiation(uint[] xValue, uint[] e, uint[] m)
        {
            uint[] x = (uint[])xValue.Clone();

            Array.Resize(ref x, m.Length);

            // mQ = -m^1 mod b
            uint mQ = Inverse(m[0]);

            // 1. temp = Mont(x, R^2 mod m), A = R mod m.
            uint[] r2 = new uint[m.Length * 2 + 1];
            r2[r2.Length - 1] = 1;

            ArithmeticKernel.Remainder(ref r2, m);

            uint[] temp = MontgomeryMultiplication(x, r2, m, mQ);

            uint[] a = new uint[m.Length + 1];
            a[a.Length - 1] = 1;

            ArithmeticKernel.Remainder(ref a, m);

            int expLen = ArithmeticKernel.Significance(e); // значащая длина exponent
            int expBits = expLen * BITS; // количество битов степени
            int expPos = expLen - 1; // позиция
            uint bitMask = BIT_MASK; // битовая маска

            // Узнаем количество значащих битов степени
            while (0 == (e[expPos] & bitMask))
            {
                bitMask >>= 1;
                expBits--;
            }

            // 2. For i from t down to 0 do the following:
            while (expBits != 0)
            {
                // 2.1 A = Mont(A, A).
                a = MontgomeryMultiplication(a, a, m, mQ);
                
                // 2.2 If e[i] = 1 then A = Mont(A, temp).
                if (0 != (e[expPos] & bitMask)) // если бит равен 1
                    a = MontgomeryMultiplication(a, temp, m, mQ);

                expBits--;
                bitMask >>= 1;

                if (0 == bitMask)
                {
                    bitMask = BIT_MASK;
                    expPos--;
                }
            }

            // 3. A Mont(A, 1).
            uint[] one = new uint[m.Length];
            one[0] = 1;

            a = MontgomeryMultiplication(a, one, m, mQ);

            ArithmeticKernel.Normalize(ref a);

            // 4. Return (A).
            return a;
        }
    }
}